
# ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
# DISCLAIMER AND GENERAL INFORMATION
#
# File: 02_map.R
# Purpose: Produces the map in section SI2
# Date: June 2025
# Data: pulled through 00_data_prep.R
#
# See 00_data_prep.R for technical disclaimer on R versions
# ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++

# Load data
source("00_data_prep.R")

# ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
# (A) Plot UK map (Figure A1) ----
# ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++

# Read different shape files
# Shapefile for local authority districts
district_shp <- st_read(dsn= "./data/Local_Authority_Districts_December_2021_GB_BUC", layer="LAD_DEC_2021_GB_BUC")  

# Shapefile for counties
county_shp <- st_read(dsn= "./data/Counties_December_2021_EN_BUC", layer="CTY_DEC_2021_EN_BUC")

# Shapefile countries
country_shp <- st_read(dsn= "./data/Countries_December_2021_UK_BUC", layer="CTRY_DEC_2021_UK_BUC")

# Convert data into simple features objects
# Districts
gbdist <- st_as_sf(district_shp, geom = sfc)
st_transform(gbdist, 4326) # Coordinate Reference System 4326
gbdist$place <- gbdist$LAD21NM # Merging variable

# Counties
gbcty <- st_as_sf(county_shp, geom = sfc) 
st_transform(gbcty, 4326) 
gbcty$place <- gbcty$CTY21NM

# Countries
gbctry <- st_as_sf(country_shp, geom = sfc) 
st_transform(gbctry)
gbctry$place <- gbctry$CTRY21NM

gbctry <- subset(gbctry, CTRY21NM != "Northern Ireland")

# Yorkshire districts
yrkdistricts <- subset(gbdist,  LAD21NM == "East Riding of Yorkshire" | 
                                LAD21NM == "Sheffield" | LAD21NM == "Doncaster" | 
                                LAD21NM == "Rotherham" | LAD21NM == "Barnsley" | 
                                LAD21NM == "Leeds" | LAD21NM == "Wakefield", 
                                select = -c(LAD21CD, LAD21NM, LAD21NMW))

# Cumbria
cumbria <- subset(gbcty, CTY21NM == "Cumbria", 
                  select = -c(CTY21CD, CTY21NM))

# Scotland and Wales
scotwales <- subset(gbctry, CTRY21NM != "England", 
                    select = -c(CTRY21CD, CTRY21NM, CTRY21NMW))

# Combine data on all target districts into single object
target <- rbind(yrkdistricts, cumbria)
target <- rbind(target, scotwales)

# Create variable classifying target areas as 'coal country', 'Scotland' or 'Wales'
target <- mutate(target, trgt = ifelse(place == "Cumbria" | 
                                       place == "Sheffield" | place == "Doncaster" | place == "Rotherham" | 
                                       place == "Barnsley" | place == "Leeds" | place == "Wakefield" | 
                                       place == "East Riding of Yorkshire", "Northern England", NA))

target <- mutate(target, trgt = ifelse(place == "Scotland", "Scotland", target$trgt))
target <- mutate(target, trgt = ifelse(place == "Wales", "Wales", target$trgt))

# Create labels for Cumbria and Yorkshire districts
county_labels <- data.frame(
  lng = c(-0.2,-3.4),
  lat = c(54,54.3),
  label = c("Yorkshire", "Cumbria"))

county_labels_sf <- st_as_sf(county_labels, coords = c("lng", "lat"), dim ="XY",
                             remove = FALSE, crs = 4326)


# Plot map
gbmap <- ggplot(data = gbctry) +
  geom_sf(data = gbctry, aes(), colour = gray(.5)) +  # specifies baseline Britain map
  geom_sf(data = target, aes(fill = trgt), color = cols[1]) +  # specifies target areas as basis for colouring
  scale_fill_manual(name="", values=c(cols[2], cols[3], cols[4]),
                    labels=c("Northern England", "Scotland", "Wales")) +
  theme(panel.background = element_rect(fill = "white"),  # sets ocean colour as white ("azure" is a good alternative)
        #legend.background = element_rect(fill=NA, linewidth=.1, linetype="solid"),  # sets legend background as transparent
        legend.position = "bottom", legend.text = element_text(size = 12),  # specifies legend position on map
        legend.key.size = unit(.5, "cm")) +  # specifies size of key boxes in legend
  geom_sf(data = county_labels_sf) +
  geom_text_repel(data = county_labels_sf, aes(x = lng, y = lat, label=label),
                  size = 4, nudge_x = c(1.4, -1.7), nudge_y = c(.9, -.5)) +  # inserts county labels and positions them
  xlab("") + ylab("") +  # removes axis labels ("lat" and "long")
  coord_sf(ylim = c(49.8, 60.5),crs= 4326)
gbmap  # displays map preview

# ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
# (B) Yorkshire inset ----
# ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++

# Add red box around the Yorkshire districts on map
gbmap2 <- gbmap +
  geom_rect(xmin = -1.9, xmax = 0.2, ymin = 54.2, ymax = 53.27, fill = NA, colour = "red", linewidth = .75) 

# Yorkshire target districts
target <- mutate(target, ykshire = ifelse(place == "East Riding of Yorkshire" | place == "Sheffield" | place == "Doncaster" | 
                                          place == "Rotherham" | place == "Barnsley" | place == "Leeds" |  place == "Wakefield", target$place, ""))

# Yorkshire districts only
yrkshire <- subset(target, ykshire != "")

# Plot map
ykmap <- ggplot(data = gbctry) +
  geom_sf(data = gbctry, color = cols[1]) +
  geom_sf(data = target, aes(fill = trgt), color = gray(.5)) + 
  scale_fill_manual(name="", values=c(cols[2], cols[3], cols[4])) +
  theme(panel.background = element_rect(fill = "white"),
        panel.border = element_rect(fill = NA), legend.position = "none") +
  xlab("") + ylab("") +
  coord_sf(xlim = c(-1.9,0.2), ylim = c(54.2,53.27), crs = 4326, expand = FALSE) +
  geom_text(data = yrkshire, aes(LONG,LAT,label = ykshire), size = 3,
            nudge_x = c(.15,0,0,.03,-.02,0,0), nudge_y = c(0,.01,.01,.02,.01,0,0), col = "black")
ykmap


# Create plot without legend for joint plotting
gbmap_wo <- gbmap2 + theme(legend.position = "none")

# Combine UK and Yorkshire maps
maps <- plot_grid(gbmap_wo, ykmap, ncol=2)

legend <- get_plot_component(gbmap, "guide-box-bottom", return_all=TRUE)

# Add legend to create full plot
plot.full <- plot_grid(maps, legend, ncol=1, rel_heights=c(1,0.1))
plot.full

# Export to PDF
ggsave(plot.full, filename="./plots/SI2_UK_map.pdf", width=10, height=8)


# ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
#                         END OF FILE
# ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++


